An intrinsic nonlinear scale governs oscillations in rapid fracture 
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When branching is suppressed, rapid cracks undergo a dynamic instability from a straight to an 
oscillatory path at a critical velocity v c . In a systematic experimental study using a wide range of 
different brittle materials, we first show how the opening profiles of straight cracks scale with the 
size £ n i of the nonlinear zone surrounding a crack's tip. We then show, for all materials tested, 
that v c is both a fixed fraction of the shear speed and, moreover, that the instability wavelength 
is proportional to l n i . These findings directly verify recent theoretical predictions and suggest that 
the nonlinear zone is not passive, but rather is closely linked to rapid crack instabilities. 



Since their discovery, a fundamental understanding of 
the origin of rapid crack instabilities [lHH has proven to 
be very elusive. The dynamics of single straight cracks 
are well described mWi by Linear Elastic Fracture Me- 
chanics (LEFM) H [lfj . This theoretical framework pre- 
dicts singular crack tip fields and describes a crack's 
dynamics as a balance between the elastic energy flux 
into the tip region and the energy dissipated at the tip. 
LEFM, however, cannot explain rapid crack instabili- 
ties and accompanying non-trivial crack patterns, with- 
out additional assumptions or physical insights about the 
near-tip region where linear elasticity breaks down. 

There have been a number of notable attempts to de- 
scribe crack instabilities by supplementing or extending 
LEFM in various ways. These include phase-field models 
(IlT - il6T |. cohesive-zone models fl7T 4l9| . models based on 
the "Principle of Local Symmetry" I20l- 2211 . energy con- 
servation bounds on crack branching [23l . |24| and models 
based on non-linear constitutive behavior near the crack 
tip [25|, .26] . Although many of these models are qualita- 
tively consistent with both experimental and numerical 
observations, decisive quantitative experiments that are 
able to differentiate between them are lacking. 

Here we focus on the oscillatory instability in rapid 
brittle fracture in which a straight crack becomes un- 
stable to sinusoidal path oscillations Q. The onset of 
these oscillations was observed at a critical velocity v c of 
about 90% of the shear wave speed c s , when the micro- 
branching instability 0, Q was suppressed. This insta- 
bility is particularly intriguing since it involves a finite 
instability wavelength at onset that is independent of ci- 
ther system geometry or loading conditions. This sug- 
gests the existence of an intrinsic scale that can not exist 
in linear elastic solutions for cracks, which are scale-free. 

Recently, a theory describing this instability was pro- 
posed [26}. This theory is based on the existence of a 
dynamic non-linear lengthscale £ n i, where linear elastic- 
ity breaks down and material nonlinearities become sig- 
nificant due to the large deformation near a crack's tip 
j27l - l30j . The basic idea behind this approach is that in 
the presence of a finite £ n i, causality implies that the 
singular LEFM fields lag behind the actual tip location 
with a delay of t& cx l n i . This led to a high- velocity os- 
cillatory instability with the following properties: (i) the 



scaled critical velocity for the onset of oscillations v c /c s 
is material independent (ii) the oscillation wavelength, 
X osc , is proportional to l n \. 

In this Letter we investigate the rapid fracture of a 
variety of different brittle gels, whose mechanical prop- 
erties vary over a wide range. We first demonstrate that 
the opening profiles of straight cracks collapse onto a sin- 
gle velocity-dependent form, when scaled by the size of 
the nonlinear elastic zone, as predicted by [29|. We then 
show that the oscillatory instability is triggered in each 
material at the same scaled value of v c /c s and, more- 
over, that the instability wavelength indeed scales with 
£ n l, confirming the theoretical predictions of [26j . 

Our experiments were performed using polyacrylamidc 
gels which are transparent, homogeneous, brittle, incom- 
pressible elastomers. The dynamics of rapid cracks in 
these neo-Hookean materials are identical to those ob- 
served in other brittle amorphous materials (e.g. glass, 
PMMA) . Due to the low elastic moduli of these soft mate- 
rials, the wave speeds and corresponding crack velocities 
are nearly 3 orders of magnitude [3~H lower than in con- 
ventional materials. This enables us to slow down the 
fracture process while obtaining detailed measurements 
of rapid cracks at unprecedented scaled velocities. 

We control the gels' physical properties by varying 
their chemical composition [32] . We varied the total 
monomer concentration (by weight) between 14.2% — 
32.4%, cross-linker concentration between 2.7% — 4.6% 
and polymer initiators in the range 0.03% — 0.06%. In 
what follows we will label each gel by its shear modulus, 
/J, (33 < fj, < 187 kPa) and fracture energy at the criti- 
cal velocity r(w c ) (24 < r(u c ) < 60 J/m 2 ). T is defined 
as the amount of energy dissipated per unit crack ex- 
tension and sample thickness. The details of these gel 
compositions are provided in [33| . Typical dimensions of 
the gels used were (x x y x z) (130 x 130 x 0.2)mm and 
(200 x 200 x 0.2)mm, where x, y and z are, respectively, 
the propagation, loading and thickness directions. 

Experiments were performed as in Q by imposing uni- 
axial tensile loading via constant displacement in the ver- 
tical (y) direction. Once a desired strain e was reached, 
a guillotine was used to initiate fracture at the sample's 
edge, midway between the vertical boundaries. To negate 
boundary effects on the crack tip prior to the onset of 



the oscillatory instability, the applied strain levels were 
selected such that the crack could reach very high ve- 
locities (~ 0.9c s ) before traversing half of the sample's 
size. For experimentally feasible system sizes, this entails 
strains e in the range 6—18%. The crack tip opening dis- 
placement (CTOD) of the moving crack was measured 
with a high speed camera focussed on an (x x y) area of 
60x9. 5-1 9mm with 1280 x 200-400 pixel resolution. Suc- 
cessive photographs were taken at between 2490/15000 
frames/s with a 2/xs exposure time. Multiple exposures 
were utilized, when needed. The micro-branching insta- 
bility was suppressed (as in Q) by setting the gel thick- 
ness to 160 — 220/im. In all experiments analyzed, no 
micro-branches occurred in the regions of interest. Post- 
fracture xy profiles were measured via an optical scanner 
with 300dpi resolution. 

Let us now consider a simple straight crack moving at 
velocity v under constant tensile loading, prior to any in- 
stability. According to LEFM, the CTOD has a parabolic 
shape whose curvature a(v) is inversely proportional to 
the instantaneous value of T(v) @. This characteris- 
tic parabolic form is indeed experimentally measured at 
points that are at a distance not too close to the crack tip. 
Sufficiently near the crack tip, regions of very high strain 
are encountered. The resulting nonlinear elastic effects 
shift the actual crack tip by a distance 6 from the apex 
of the parabolic form defined by LEFM [27|. The dissi- 
pative zone adjacent to the tip is also contained within S. 
In gels, the dissipative zone is significantly smaller than 
the size of the non- linear elastic deformation zone [28j |. 

The strain levels imposed in our measurements suggest 
that the CTOD predicted by LEFM should be calculated 
with respect to the background strain e. To this end, we 
consider the energy functional describing our incompress- 
ible gels under plane stress conditions |34j 
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where F — Wip is the deformation gradient and x' = ip(x) 
is a mapping between a reference (undeformed) configu- 
ration x and a deformed one x' . For our uniaxial loading 
we have tp x = (1 + e)~ 1 ^ 2 x + u x and (p y = (1 + ()y + u y , 
where u is the displacement field due to the presence of a 
crack. Using the stress measure s = 8fU(F), the momen- 
tum balance equation reads V • s = pduf, where p is the 
mass density. The traction-free boundary conditions on 
the crack faces take the form s xy (r, 6 = ±w) = s yy (r, 6 = 
±7r) = 0, where (r, 6) is a polar coordinate system moving 
with the crack tip and = is the propagation direction. 
Linearizing these equations with respect to u and solv- 
ing the resulting equations numerically near the tip of a 
crack mo ving at a steady velocity v, the solution takes 
the form 13311 
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FIG. 1. (Color online) Comparison of the fracture energy, 
T(d), calculated from (squares) CTOD measurements using 
the LEFM solution and (circles) using the extension of 
LEFM for finite strain given by Eq. ©. The gel used has 
(j, — HOkPa. (inset) Measurements of T(v) obtained by ve- 
locity profiles in a strip geometry with a more compliant 
(^i = 36kPa) gel 0. 



where Kj is the stress intensity factor, T is a traction 
vector known as the T-stress and $1 is a universal angular 
function HH . In the limit e -> the s tandard LEFM 
solution is recovered Q. Equation (|A15|) can be used to 
relate the parabolic crack tip curvature a(v,e) and the 
fracture energy T(v), yielding 
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where T x (e) and A(v; e) are given in [331 ] . 

The deviations of Eq. §3§ from the LEFM predictions 
are very small for low strains (e < 6%), justifying their 
omission in previous studies [27H30f . However, as e is in- 
creased, significant corrections to the LEFM predictions 
appear for high crack velocities. In Fig. [TJ we compare 
values of T(v) derived from measured CTOD's using both 
Eq. ([3]) and the analogous LEFM relation. While, at low 
velocities, the differences arc insignificant, for v > 0.85c s 
the T(v) curves diverge significantly. The LEFM relation 
(blue-shading in Fig. [1]) yields both a large systematic 
variation of T(v) for different e values and a substantial 
decrease of T with v. The former contradicts the expec- 
tation that T(v) is a material-dependent function whose 
value should not reflect the background strain. Moreover, 
the sharp decrease in T(v) at high velocities is surprising. 
One would expect a nearly constant value of T(v) in this 
narrow range of velocities, as obtained in independent 
measurements of T(v) (using an "infinite strip" geome- 
try Q) for a similar material (Fig. [TJ - inset). The use 
of Eq. ^ both eliminates the apparent dependence of T 
with e and indeed reveals the slowly increasing T with v, 
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consistent with the direct measurements presented in the 
inset. We were able to significantly decrease the 2 — 4% 
experimental uncertainty in v/c s , by varying c s (within 
experimental uncertainty) to minimize the variance of 
the mean value of T(v) over the range 0.85 <v/c s < 0.9. 
r(w), obtained by this procedure, is presented within the 
red-shaded data in Fig. [T] The collapse of the data to- 
gether with the resulting slow increase of T(v) with v as 
expected from Q, justifies this procedure, which is used 
to determine T(v) in what follows. 

For each v/c s , Fig. [5^, demonstrates that scaling 
lengths by T/fj, collapses the CTOD's of different ma- 
terials to a single function. LEFM predicts that this 
should occur for the parabolic CTOD's away from the 
crack tip. It is not obvious, however, that data collapse 
should occur in the near-tip region defined by S, as this 
is a wholly independent regime. Data collapse with T//J, 
in the weakly nonlinear regime (i.e. cubic expansion of 
U in the metric strain measure E = ^(F T F—I) (HI) was 
predicted for nco-Hookean materials [28l . I29I ] , where sec- 
ond order elastic coefficients are order fi. (In analogous 
scaling for other materials these coefficients may signifi- 
cantly differ from fi [33j].) A perfect data collapse would 
indicate that this is the only significant scale in the sys- 
tem. High- resolution measurements of S(v), presented in 
Fig. [2b for 5 different materials, provide a stringent test 
of this scaling. While the widely spread raw data (Fig. 
[5}>top) indeed undergo an approximate collapse when 
scaled by V/fi, the imperfect collapse for small values of 
scaled S indicates that an additional, much smaller, scale 
exists. We surmise that this additional scale could be 
related to either the strongly nonlinear elastic region or 
the dissipation zone |28| . 

We now turn to the oscillatory instability. As shown 
in Fig. [3j the wavelength of the first oscillation, X osc , is 
strongly material-dependent, varying by over a factor of 
2.5 in different materials (Fig. [Sb). In each material there 
is a well-defined velocity v c for the onset of the instability. 
As predicted by [111, Fig. shows that v c , when scaled 
by c s , has the nearly constant value of w c = 0.9c s , in each 
of the 6 materials studied. 

What is the origin of the instability wavelength? Fig. 
[3] confirms that X osc is not related to details of the experi- 
mental system. In experiments with identical conditions, 
X osc varied widely with the material used. In (2(|, X osc 
was predicted to be proportional to the size of the nonlin- 
ear zone £ n i. Here we use 5{v = v c ) to estimate £„i at the 
critical velocity v c for different materials. The obvious 
advantage of doing this is that S(v) is directly measurable, 
and hence the theoretical prediction of [26| can be recast 
as a relation between two directly measurable quantities, 
A osc and 5{v = v c ). In Fig. [4] we plot A osc vs. S(v = v c ) for 
the 6 materials used. We indeed find that S is directly 
proportional to A osc , as predicted in (26j. Moreover, the 
constant of proportionality between A osc an d 5( v = v c ) in 
Fig. [4] is consistent with the analysis of [26l |33|. 

We note that the weakly nonlinear estimate of £ n i cx 
T/n [2!| is also linearly related to A osc . In contrast to 
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FIG. 2. (Color online) (a) Top: Measurements of the CTOD 
for 5 different materials (legend below) at £ ~ 0.08 and 
u = 0.86c s . (x',y') are the coordinates in the laboratory (de- 
formed) frame, to be distinguished from the reference (un- 
deformed) frame (x,y). Bottom: When scaled by T(v)/fj,, 
these curves collapse to a single function. Far from the tip 
at (x' ,y') = (0,0) the CTOD is parabolic (dashed line), but 
strongly deviates from this form at a scale S, defined as the 
distance from the apex of these parabola to the crack tip. (b) 
Top: 8(v/cs) for the 5 different materials in (a). Bottom: 
approximate collapse of the 8(v/c s ) curves when scaled by 

r»//i . 




FIG. 3. (Color online) (a) Typical photographs of the xy 
profiles of fracture surfaces at the onset of the oscillatory in- 
stability; from top to bottom: /i = 36kPa, 143kPa, 168kPa. 
(b) The oscillation wavelength, A osc changes significantly with 
the material, as characterized by [i. (c) The scaled critical ve- 
locity, v c — 0.9c s , is constant. v c is defined as the maximal 
velocity prior to the instability onset in each material. Symbol 
colors correspond to the legend of Fig. [2] 
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FIG. 4. (Color online) Comparison between the non-linear 
length-scale S(v = v c ) and the oscillation wavelength X oac . 
Note that the different combinations of /j,, T, and e are used 
to produce ~ 15 independent measurements. The dashed line 
is a guide to the eye. 



Fig. |4l however, this linear relation involves an offset cor- 
responding to ~ 100 — 300/im. This scale is also apparent 
in the imperfect data collapse in Fig. ^jp, suggesting that 
d includes length-scales such as the strongly nonlinear 
contributions to the nonlinear elastic zone and/or the 



scale of the "dissipative zone" at the crack tip which are 
beyond the perturbative estimate of Eq. (JTJ used in (2t| . 

In conclusion, our results conclusively demonstrate 
that the oscillatory instability of fast brittle cracks in- 
deed involves an intrinsic scale that is governed, in a 
large part, by the nonlinear elastic zone surrounding the 
crack tip. The size of this zone quantitatively agrees 
with the predictions of [26j. These results indicate that 
the nonlinear (and dissipative) zones surrounding the tip 
of a moving crack are not "passive" objects that are sim- 
ply "dragged along" by the crack tip. Instead, as sug- 
gested by [1, HH, Hi| , this region may play an active role 
in destabilizing crack motion. The demonstration of this 
presented in this work is, therefore, an important step 
in obtaining a fundamental understanding of the origin 
of instabilities in dynamic fracture. These ideas are as 
general as the singular behavior that occurs at the tip of 
a moving crack. It is therefore conceivable that dynam- 
ics of the near-tip zone could play an important role in 
unraveling the physical mechanism driving other insta- 
bilities of rapid cracks [TMBI- 
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I. MATERIALS 



Our experiments were carried out on cross- linked polyacrylamide gels. Their elastic properties are determined 
by the molecular weight of the monomer (acrylamide) chains and by the concentration of cross-linking molecules 
(bis-acrylamidc). The gels were prepared by adding to an acrylamidc/bis-acrylamide solution (with varying concen- 
trations) two initiators: ammonium persulfate (APS) and tetramethyl ethylene diamine (TEMED). We used a fixed 
0.63% (weight/ Volume) concentration of the APS while changing the TEMED concentration, thus controlling the 
length of polymers chains and their respective molecular weight. The following table provides a detailed description 
of the gels' composition and their measured shear moduli: 





Total Monomer (%) 


Cross-Linker/Monomer (%) 


Total Initiator (TEMED) (100 x %) 


Shear Modulus, fi (kPa) 


1 


14.1 


2.7 


6.2 


36 ±0.5 


2 


27 


4.6 


6.2 


168 ±3 


3 


22.7 


4.2 


6.2 


110 ±6 


4 


31.4 


3.2 


6.2 


143 ± 14 


5 


14.1 


2.7 


3.1 


34 ± 1.5 


6 


27 


4.6 


3.1 


187 ± 17 



II. DERIVATION OF EQS. (2)-(3) IN THE MAIN TEXT 



The polyacrylamide gels in our experiments are incompressible neo-Hookean materials deformed under plane stress 
conditions. They are described by the following energy functional [Sll | 

U(F) = I (FijFij + [det(F)}- 2 3) , (A4) 

where F is the 2D deformation gradient F = Wip and tp is the deformation. 

As explained in the main text, experimental limitations on the system sizes we could use force us to impose prestrains 
e along the y-direction such that the deformation (p(x) takes the form 

fx = (1 + e)~ 1/2 x + u x , ip y = (1 + e)y + u y . (A5) 
The transverse stretch A is given by 

A = [det(F)]- 1 = J— . (A6) 

d X ip x dyiPy ~ d X (PyC)ytp x 

Note that when u = 0, we have A = d x tp x = (1 + e)~ 1//2 , as expected for an isotropic incompressible material under 
homogeneous uniaxial loading. The first Piola-Kirchhoff tensor s, derived from s-dpU(F), reads 

Sij = [i (djtpi - A 3 Ei k e jl di(pk) , (A7) 

where e,-j is the 2D alternator, i.e. e xx =£ yy = 0, e xy = — e yx — 1. The equations of motion V ■ s = pip can be written 
explicitly as 

fiV 2 if X + [l [dy(A 3 )d X <Py - d X (A 3 )dy<Py] = fi f X , (A8) 

fiV 2 (p y + n [d x (A 3 )d y ip x - d y (A 3 )d x <p x ] = pcp y , (A9) 



G 



where p is the reference mass density. The traction- free boundary conditions on the crack faces s xy (r,8 = ±7r) 
s yy (r, 9 — ±tt) =0 can be written explicitly as 



s xy (r, 9 = ±7r) = fi(d y ip x + A 3 d x <p y )\e= 



±7T 



0. 



Syy(r, I 



±7T) 



P(dytp y ~ A 3 d x tp x )\g-_ 



(A10) 



where the polar coordinate system (r, 0) is moving with the crack tip such that 9 = is the crack tip propagation 
direction. 

Focusing on steady state conditions, d t = —vd x , we can now linearize Eqs. (|A8[) . (|A9|) and (|A10[) with respect to u. 
The results read 



49««i + 3(1 + e) 3/2 d xy u y + d yy u x - v 2 d x 



, 



d xx u y + 3(l + e)- 3/2 d xy u x + 1 + 3 (1 + e) 



with the following boundary conditions at 



v 

±7T 



OyyUy ~ ' 8 XX Uy = , 



l + 3(l + e) 



-3 



d y U X + (1 + e)- 3 ' 2 d X Uy = , 

d y u y + 2(l + e)- 3 ^d x u x =0 , 



(All) 
(A12) 

(A13) 
(A14) 



where v here is normalized by c s = \fp~Jp. In the limit e — > 0, we recover the equations of isotropic LEFM. The 
two-term asymptotic expansion near the crack tip appears in Eq. (2) in the main text and is copied here 



u x (r,9;v,e) = ^fej^ O^fl; v, e) + 



T x (e)rcos6> 



Uy(r,6;v,e) 
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The function O(0; v, e) is calculated using the following half-integer Fourier expansion 

"(2n-l) 



tt x (9; v, e) = a„(u, e) cos 



17^(0; i>, e) = &„(d, e) sin 



(2n- 1)0 



(A15) 



(A16) 



{a n , b n } are calculated by solving a set of linear algebraic equations for each v and e, where n is selected such that the 
required accuracy is obtained. Note that since this is a linear set of equations, {a n ,b n } can be determined up to an 
overall multiplicative factor, which as usual is quantified by the stress intensity factor Kj that cannot be calculated 
by the asymptotic analysis (but rather from the global crack problem). The standard LEFM result is recovered for 
e — > 0. For the sub-leading term we obtain 



T x (e) 



T( £ )[l + 3(l + 6)- 
12 



Ty(e) 



T(e) (l + e)- 3/2 



6 



(A17) 



where T(e) cannot be determined by the asymptotic analysis. Note that, again, the standard LEFM result is recovered 
for e — > 0. In the absence of a global nonlinear solution of our problem, we use the latter to obtain a sensible estimate 
of T x (e), which appears in Eq. (3) in the main text. It is known that for a tensile crack in a large system under 
remote tensile stress er°°, T x (e)/fi = —tj°°/E, where E is the Young's modulus |S2j ] . Since the background strains in 
our experiments were of the order of 0.1, we use the linear approximation cr°°~ei? to obtain T x (e)/n~—e. We used 
this estimate in our analysis. 

To obtain the crack tip opening displacement (CTOD) we focus on 9 = ir, which leads to 



n y (ir;v,e) 



<p x (r,ir) = -r\(l + e) 1/2 - el , (p y (r,ir) = ^^Z- 

L J Afly/ZTT 

Eliminating r between the last two relations, we obtain a parabolic form (f x (r, w) = — a ip 2 (r, ir) with 

32717? 



(1 



-1/2 



nl(7r;v,e)K 2 



The final step in deriving Eq. (3) in the main text would be to consider the J-intcgral [S5 

1 



G = 



U + -pv d x Uid x Ui ) n x - SijrijdxUi 



dC 



(A18) 



(A19) 



(A20) 
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where C is any path surrounding the crack tip. We calculated the integral numerically using Eqs. (IA15[) . which allows 
us to define A(v; e) from the relation 



A( v;e)K 2 ( v,e) 
32tt^ 2 



G ^ICT 1 - (A21) 



Using energy balance, T(v) = G, together with Eq. (| A19|) . we arrive at Eq. (3) in the main text. 



III. ADDITIONAL COMMENTS 



• The constant of proportionality between A osc and S(v c ) 

The constant of proportionality between X osc and S(v c ) in Fig. 4 in the main text is about 10. To understand 



this, consider Eq. (14) in |S4 



P -Mv c ), (A22) 



where c n i is the typical wave speed within the nonlinear zone and j3 is a dimensionless number quantifying the 
typical time it takes mechanical information to propagate across the nonlinear zone (and is order unity). The 
real part of the dimensionless complex oscillation frequency was estimated to be 5i(<D c ) ~ 1.5. Furthermore, 
/3~l/2 and c„;~c s was used, leading to A osc ~ 1.6£ n i(v c ) [S4| . Finally, £ n i(v c ) was estimated using a weakly 
nonlinear calculation ahead of the crack tip [S|| to be of the order of a few mm. The relation A osc ~ 10 8(v c ) may 
emerge from two observations: (i) S(v c ) is measured behind the crack tip, that may be somewhat smaller than the 
estimate ahead of the tip, leading to £ n i (v c ) = a 8(v c ), with a > 1 (ii) /3 may be more faithfully estimated as (3 = 2, 
which corresponds to a round trip across the nonlinear zone, increasing the prefactor 1.6 in A osc ~ 1.6£ n i(v c ) to 
about 6. 

• The nonlinear scale T/fi 

In the main text the weakly nonlinear scale T/fx is used, e.g. Fig. 2. We would like to stress that in the 

most general case a weakly nonlinear elastic lengthscale depends on a set of second order elastic coefficients, 
in addition to the linear elastic coefficient fi and A (the first Lame constant). To see this, write the following 
expansion of the energy functional U in 2D 

U 2D (E) = -{trEf +[itrE 2 + -(trEf + —trE 3 , (A23) 
2 3 3 

where the metric (Green-Lagrange) strain measure E= \{F T F — I) is understood as a 2D tensor and {l,m} 
are two of the Murnaghan coefficients [S6|. In 3D there is an additional second order coefficient n. Eq. (|A23|) 
reduces to Eq. (1) in the main text when A = 2fi, I = — 4/i, m = — 4/x. For other materials I and m may be 
larger. For example, for Polystyrene (a glassy polymer) 14/i [S7| |. Since the weakly nonlinear estimate for 
£ n i involves (I //i) 2 [S8I |. T//J, (which is of course dimensionally correct) may include widely varying prefactors 
for different materials. 
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